load required packages

library(qiime2R)
library(tidyverse)
library(ggplot2)
library(lubridate)
library(xts)
library(ggpubr)
library(dplyr)
library(data.table)
library(vegan)
#library(reshape2)
#install.packages("janitor")
#library(janitor)
library(viridis)

show versions:

sessionInfo()
## R version 4.3.2 (2023-10-31)
## Platform: x86_64-apple-darwin20 (64-bit)
## Running under: macOS Monterey 12.2.1
## 
## Matrix products: default
## BLAS:   /Library/Frameworks/R.framework/Versions/4.3-x86_64/Resources/lib/libRblas.0.dylib 
## LAPACK: /Library/Frameworks/R.framework/Versions/4.3-x86_64/Resources/lib/libRlapack.dylib;  LAPACK version 3.11.0
## 
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
## 
## time zone: Europe/London
## tzcode source: internal
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] viridis_0.6.4     viridisLite_0.4.2 vegan_2.6-4       lattice_0.21-9   
##  [5] permute_0.9-7     data.table_1.15.0 ggpubr_0.6.0      xts_0.13.2       
##  [9] zoo_1.8-12        lubridate_1.9.3   forcats_1.0.0     stringr_1.5.1    
## [13] dplyr_1.1.4       purrr_1.0.2       readr_2.1.4       tidyr_1.3.1      
## [17] tibble_3.2.1      ggplot2_3.4.4     tidyverse_2.0.0   qiime2R_0.99.6   
## 
## loaded via a namespace (and not attached):
##  [1] bitops_1.0-7            gridExtra_2.3           rlang_1.1.3            
##  [4] magrittr_2.0.3          ade4_1.7-22             compiler_4.3.2         
##  [7] mgcv_1.9-0              vctrs_0.6.5             reshape2_1.4.4         
## [10] pkgconfig_2.0.3         crayon_1.5.2            fastmap_1.1.1          
## [13] backports_1.4.1         XVector_0.42.0          utf8_1.2.4             
## [16] rmarkdown_2.25          tzdb_0.4.0              xfun_0.42              
## [19] zlibbioc_1.48.0         cachem_1.0.8            GenomeInfoDb_1.38.4    
## [22] jsonlite_1.8.8          biomformat_1.30.0       rhdf5filters_1.14.1    
## [25] Rhdf5lib_1.24.1         broom_1.0.5             parallel_4.3.2         
## [28] cluster_2.1.4           R6_2.5.1                bslib_0.6.1            
## [31] stringi_1.8.3           car_3.1-2               zCompositions_1.5      
## [34] rpart_4.1.21            jquerylib_0.1.4         Rcpp_1.0.12            
## [37] iterators_1.0.14        knitr_1.45              base64enc_0.1-3        
## [40] IRanges_2.36.0          Matrix_1.6-1.1          splines_4.3.2          
## [43] nnet_7.3-19             igraph_1.6.0            timechange_0.2.0       
## [46] tidyselect_1.2.0        abind_1.4-5             rstudioapi_0.15.0      
## [49] yaml_2.3.8              codetools_0.2-19        plyr_1.8.9             
## [52] Biobase_2.62.0          withr_3.0.0             evaluate_0.23          
## [55] foreign_0.8-85          survival_3.5-7          Biostrings_2.70.1      
## [58] pillar_1.9.0            phyloseq_1.46.0         carData_3.0-5          
## [61] checkmate_2.3.1         DT_0.31                 foreach_1.5.2          
## [64] stats4_4.3.2            NADA_1.6-1.1            generics_0.1.3         
## [67] RCurl_1.98-1.13         truncnorm_1.0-9         S4Vectors_0.40.2       
## [70] hms_1.1.3               munsell_0.5.0           scales_1.3.0           
## [73] glue_1.7.0              Hmisc_5.1-1             tools_4.3.2            
## [76] ggsignif_0.6.4          rhdf5_2.46.1            grid_4.3.2             
## [79] ape_5.7-1               colorspace_2.1-0        nlme_3.1-163           
## [82] GenomeInfoDbData_1.2.11 htmlTable_2.4.2         Formula_1.2-5          
## [85] cli_3.6.2               fansi_1.0.6             gtable_0.3.4           
## [88] rstatix_0.7.2           sass_0.4.8              digest_0.6.34          
## [91] BiocGenerics_0.48.1     htmlwidgets_1.6.4       htmltools_0.5.7        
## [94] multtest_2.58.0         lifecycle_1.0.4         MASS_7.3-60

Eukaryotes

18S metabarcoding data analysis ## Beta Diversity

setwd("/Users/RMS1U18/Library/CloudStorage/OneDrive-UniversityofSouthampton/PhD/Soton/Bioinformatics/18S_dada2redo/")

SVs <- read_qza("rarefied_table.qza")$data
metadata <- read_q2metadata("metadata-rarefied.txt")
taxonomy <- read_qza("taxonomy.qza")$data %>% parse_taxonomy()
taxonomy$Kingdom <- gsub("d__", "",taxonomy$Kingdom)

metadata <- data.frame(metadata[,-1], row.names=metadata[,1])

taxasums<-summarize_taxa(SVs, taxonomy)$Species

# now I'm going to create betadiversity plots and look at the what environmental
# variables may be influencing beta diversity
mtaxasums <- t(as.matrix(taxasums))
set.seed(123)
nmds <- metaMDS(mtaxasums, distance = "bray")
## Square root transformation
## Wisconsin double standardization
## Run 0 stress 0.1898539 
## Run 1 stress 0.1947476 
## Run 2 stress 0.1918644 
## Run 3 stress 0.1998882 
## Run 4 stress 0.194419 
## Run 5 stress 0.1952102 
## Run 6 stress 0.1916873 
## Run 7 stress 0.200839 
## Run 8 stress 0.2079919 
## Run 9 stress 0.1994643 
## Run 10 stress 0.3982293 
## Run 11 stress 0.1940355 
## Run 12 stress 0.1943068 
## Run 13 stress 0.2128468 
## Run 14 stress 0.18729 
## ... New best solution
## ... Procrustes: rmse 0.05572826  max resid 0.1936431 
## Run 15 stress 0.2062908 
## Run 16 stress 0.1991352 
## Run 17 stress 0.2159338 
## Run 18 stress 0.189858 
## Run 19 stress 0.1876837 
## ... Procrustes: rmse 0.02183372  max resid 0.09050004 
## Run 20 stress 0.1940356 
## *** Best solution was not repeated -- monoMDS stopping criteria:
##     20: stress ratio > sratmax
en <- envfit(nmds, metadata, permutations = 999, na.rm = TRUE)

#extract NMDS scores (x and y coordinates) for sites from newer versions of vegan package
data.scores <- as.data.frame(scores(nmds)$sites)
data.score <- merge(data.scores, metadata, by=0) # adding metadata
#extracting vectors
en_coord_cont = as.data.frame(scores(en, "vectors")) * ordiArrowMul(en)
en_coord_cat = as.data.frame(scores(en, "factors")) * ordiArrowMul(en)

Mantel Tests

# To work out which environmental vectors to include I'm going to do some mantel tests
# first create the distance matrices then do the mantel test
dist.abundance <- vegdist(t(taxasums), method = "bray")
dist.time <- dist(metadata$Time.Since.Start, method = "euclidean")
abund_time = mantel(dist.abundance, dist.time, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_time
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.time, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.1511 
##       Significance: 0.005 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0700 0.0931 0.1112 0.1332 
## Permutation: free
## Number of permutations: 9999
dist.day.night <- dist(metadata$nDay_Night, method = "euclidean")
abund_day.night = mantel(dist.abundance, dist.day.night, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_day.night
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.day.night, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: -0.03779 
##       Significance: 0.7517 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0688 0.0903 0.1081 0.1353 
## Permutation: free
## Number of permutations: 9999
dist.moon <- dist(metadata$nMoon_Moonset, method = "euclidean")
abund_moon = mantel(dist.abundance, dist.moon, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_moon
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.moon, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: -0.09146 
##       Significance: 0.9681 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0727 0.0956 0.1166 0.1413 
## Permutation: free
## Number of permutations: 9999
dist.tide <- dist(metadata$UKHO_Tide_Height, method = "euclidean")
abund_tide = mantel(dist.abundance, dist.tide, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_tide
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.tide, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: -0.01736 
##       Significance: 0.6145 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0638 0.0853 0.1053 0.1244 
## Permutation: free
## Number of permutations: 9999
dist.caltemp <- dist(metadata$Calista_Temp_oC, method = "euclidean")
abund_caltemp = mantel(dist.abundance, dist.caltemp, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_caltemp
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.caltemp, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.2309 
##       Significance: 0.0026 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0976 0.1295 0.1567 0.1846 
## Permutation: free
## Number of permutations: 9999
dist.airtemp <- dist(metadata$SotonMet_Air_Temp_oC, method = "euclidean")
abund_airtemp = mantel(dist.abundance, dist.airtemp, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_airtemp
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.airtemp, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.02534 
##       Significance: 0.3255 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0858 0.1127 0.1373 0.1621 
## Permutation: free
## Number of permutations: 9999
dist.wind <- dist(metadata$SotonMet_Wind_Speed_knots, method = "euclidean")
abund_wind = mantel(dist.abundance, dist.wind, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_wind
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.wind, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.132 
##       Significance: 0.0801 
## 
## Upper quantiles of permutations (null model):
##   90%   95% 97.5%   99% 
## 0.118 0.156 0.187 0.224 
## Permutation: free
## Number of permutations: 9999
dist.winddirection <- dist(metadata$SotonMet_Wind_Direction_degrees, method = "euclidean")
abund_winddirection = mantel(dist.abundance, dist.winddirection, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_winddirection
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.winddirection, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.101 
##       Significance: 0.1405 
## 
## Upper quantiles of permutations (null model):
##   90%   95% 97.5%   99% 
## 0.121 0.158 0.193 0.235 
## Permutation: free
## Number of permutations: 9999
dist.gust <- dist(metadata$SotonMet_Max_Gust_knots, method = "euclidean")
abund_gust = mantel(dist.abundance, dist.gust, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_gust
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.gust, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.07845 
##       Significance: 0.1727 
## 
## Upper quantiles of permutations (null model):
##   90%   95% 97.5%   99% 
## 0.109 0.146 0.179 0.214 
## Permutation: free
## Number of permutations: 9999
dist.pressure <- dist(metadata$SotonMet_Barometric_Pressure_millibars, method = "euclidean")
abund_pressure = mantel(dist.abundance, dist.pressure, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_pressure
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.pressure, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.1275 
##       Significance: 0.0933 
## 
## Upper quantiles of permutations (null model):
##   90%   95% 97.5%   99% 
## 0.123 0.163 0.199 0.239 
## Permutation: free
## Number of permutations: 9999
dist.UoStemp <- dist(metadata$UoS_Temp_oC, method = "euclidean")
abund_UoStemp = mantel(dist.abundance, dist.UoStemp, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_UoStemp
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.UoStemp, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.1409 
##       Significance: 0.0455 
## 
## Upper quantiles of permutations (null model):
##   90%   95% 97.5%   99% 
## 0.104 0.136 0.167 0.195 
## Permutation: free
## Number of permutations: 9999

nMDS

#removing non significant variables and time since start
drops <- c("Conc_NOC_ngµl",
           "Conc_ESS_ngµl",
           "Conc_PCR2",
           "Conc_diff_ngµl",
           "nDay_Night",
           "nMoon_Moonset",
           "SampleNo",
           "UKHO_Tide_Height",
           "SotonMet_Air_Temp_oC",
           "SotonMet_Wind_Speed_knots",
           "SotonMet_Wind_Direction_degrees",
           "SotonMet_Max_Gust_knots",
           "SotonMet_Barometric_Pressure_millibars",
           "UoS_Temp_oC", 
           "Time.Since.Start")
en_coord_cont <- en_coord_cont[!(row.names(en_coord_cont) %in% drops),]

nMDS18S = ggplot(data = data.score, aes(x = NMDS1, y = NMDS2)) + 
  geom_point(data = data.score, aes(color=`Time.Since.Start`/60), size = 3, alpha = 0.8) + 
  scale_fill_viridis(limits = c(0, 94))  +
  scale_colour_viridis(limits = c(0, 94), 
                        breaks = c(0, 24, 48, 72),
                        labels = c(0, 24, 48, 72)) +
  geom_segment(aes(x = 0, y = 0, xend = NMDS1, yend = NMDS2), 
               data = en_coord_cont, linewidth =1, alpha = 0.5, colour = "grey30") +
  geom_text(data = en_coord_cont, aes(x = NMDS1, y = NMDS2), colour = "grey30", 
            fontface = "bold", label = "Temp", parse = TRUE, position = position_nudge(x = -0.05)) + 
  theme(axis.title = element_text(size = 10, face = "bold", colour = "grey30"), 
        panel.background = element_blank(), panel.border = element_rect(fill = NA, colour = "grey30"), 
        axis.ticks = element_blank(), axis.text = element_blank(), legend.key = element_blank(), 
        legend.title = element_text(size = 10, face = "bold", colour = "grey30"), 
        legend.text = element_text(size = 9, colour = "grey30")) +
  theme_q2r() +
  labs(title="Eukaryotes", colour = "Hours Since Start")
nMDS18S

## Alpha Diversity Time Series

shannon<-read_qza("/Users/RMS1U18/Library/CloudStorage/OneDrive-UniversityofSouthampton/PhD/Soton/Bioinformatics/18S_dada2redo/core-metrics-results/shannon_vector.qza")
shannon<-shannon$data %>% rownames_to_column("SampleID") # this moves the sample names to a new column that matches the metadata and allows them to be merged

observedASVs<-read_qza("/Users/RMS1U18/Library/CloudStorage/OneDrive-UniversityofSouthampton/PhD/Soton/Bioinformatics/18S_dada2redo/core-metrics-results/observed_features_vector.qza")
observedASVs<-observedASVs$data %>% rownames_to_column("SampleID")

evenness<-read_qza("/Users/RMS1U18/Library/CloudStorage/OneDrive-UniversityofSouthampton/PhD/Soton/Bioinformatics/18S_dada2redo/core-metrics-results/evenness_vector.qza")
evenness<-evenness$data %>% rownames_to_column("SampleID")

faithPD<-read_qza("/Users/RMS1U18/Library/CloudStorage/OneDrive-UniversityofSouthampton/PhD/Soton/Bioinformatics/18S_dada2redo/core-metrics-results/faith_pd_vector.qza")
faithPD<-faithPD$data %>% rownames_to_column("SampleID")

metadata<- metadata %>% rownames_to_column("SampleID") # I think this just assigned 1,2,3,4
# instead I need to change the Sample column to SampleID
## first remove the new sample ID column
# metadata = subset(metadata, select = -SampleID)
# then rename Sample to SampleID
# names(metadata)[names(metadata) == 'sample'] <- 'SampleID'

#assign alpha diversity metrics to the metadata
metadata<-
  metadata %>% 
  left_join(shannon) %>%
  left_join(observedASVs) %>%
  left_join(evenness) %>%
  left_join(faithPD)
## Joining with `by = join_by(SampleID)`
## Joining with `by = join_by(SampleID)`
## Joining with `by = join_by(SampleID)`
## Joining with `by = join_by(SampleID)`
metadata$`Start_Date_Time`<-dmy_hm(metadata$`Start_Date_Time`) # getting the dates into ISO format

ShannonPlot <- ggplot(metadata, aes(x = `Start_Date_Time`, y=`shannon_entropy`)) +
  geom_line() +
  geom_point() +
  ylim(2,7.5) +
  theme_q2r() +
  xlab("Time") +
  ylab("Diversity")
ShannonPlot 

EvennessPlot <- ggplot(metadata, aes(x = `Start_Date_Time`, y=`pielou_evenness`)) +
  geom_line() +
  geom_point() +
  ylim(0.25,0.85) +
  theme_q2r() +
  xlab("Time") +
  ylab("Evenness")
EvennessPlot 

ObservedPlot <- ggplot(metadata, aes(x = `Start_Date_Time`, y=`observed_features`)) +
  geom_line() +
  geom_point() +
  theme_q2r() +
  xlab("Time") +
  ylab("ASVs") +
  labs(title = "Eukaryotes")
ObservedPlot

FaithPlot <- ggplot(metadata, aes(x = `Start_Date_Time`, y=`faith_pd`)) +
  geom_line() +
  geom_point() +
  theme_q2r() +
  xlab("Time") +
  ylab("Faith PD")
FaithPlot 

Alpha18STime<- ggarrange(ObservedPlot, ShannonPlot, EvennessPlot, FaithPlot,
          ncol = 1, 
          align = "hv")
Alpha18STime

## Stacked Taxa Time Series ### Species

setwd("/Users/RMS1U18/Library/CloudStorage/OneDrive-UniversityofSouthampton/PhD/Soton/Bioinformatics/18S_dada2redo/")
metadata <- read_q2metadata("metadata-rarefied.txt")
metadata <- data.frame(metadata[,-1], row.names=metadata[,1])
taxasums <- cbind(taxasums, total = rowSums(taxasums)) #add total col for each species
taxasums <- taxasums[order(-taxasums$total),] #order by total column
write.csv(taxasums, "SpeciesList.csv")

# seperate the top 10 taxa
top10species <- taxasums %>%
  arrange(desc(total)) %>%
  slice(1:10)
#remove the top 10 species 
remainderSpecies <- taxasums[-(1:10),]
remainderSpecies <- rbind(remainderSpecies, Remainder = colSums(remainderSpecies)) #add total row for each sample
# rejoin the tables
taxasums10 <- rbind(top10species, remainderSpecies)
# take the top 11 which will now include the remainder total
taxasums10 <- taxasums10 %>%
  arrange(desc(total)) %>%
  slice(1:11)
#remove the total column
taxasums10$total <- NULL
#rownames to first column
taxasums10<- tibble::rownames_to_column(taxasums10)
taxasums10<- data.frame(t(taxasums10)) #transpose
names(taxasums10)<- taxasums10[1,] # make the taxa the rownames
taxasums10<- taxasums10[-1,] # remove the first duplicate first row

# add Date_Time to taxasums10
t.metadata <- metadata[,c("sample","Start_Date_Time")]
t.metadata$Start_Date_Time <- as.character(t.metadata$Start_Date_Time) #convert posixct back to character



#t.metadata <- row_to_names(t.metadata, row_number = 1, remove_row = TRUE) #make the sample names col names
taxasums101 <- cbind(t.metadata, taxasums10) #combine taxa and times
taxasums101$Start_Date_Time<- dmy_hm(taxasums101$Start_Date_Time)
write.csv(taxasums101, "species10.csv")

taxasums101$sample<- NULL # remove sample column

species<- taxasums101 %>%
  pivot_longer(!Start_Date_Time,, names_to = "Taxa", values_to = "Sequence count")

species$`Sequence count`<- as.numeric(species$`Sequence count`)

#plot stack
species18SStack <- species %>%
  ggplot(aes(fill= Taxa,
             x=Start_Date_Time,
             y=`Sequence count`,
             colour = Taxa)) +
  geom_area(position="stack", 
            stat="identity") +
  geom_line(position = 'stack', linewidth = 1, colour = "black" ) +
  scale_colour_viridis(discrete = TRUE, direction = -1) +
  scale_fill_viridis(discrete = TRUE, direction = -1) +
  theme_q2r() +
  theme(legend.position = "bottom", legend.box = "vertical", legend.justification = "left", legend.box.just = "left", legend.direction = "vertical") +
  labs(title = "Eukaryotes", x="Time")
  

species18SStack

### Genus

Prokaryotes

Beta Diversity

setwd("/Users/RMS1U18/Library/CloudStorage/OneDrive-UniversityofSouthampton/PhD/Soton/Bioinformatics/16S/")

SVs <- read_qza("rarefied_table.qza")$data
metadata <- read_q2metadata("metadata-no-blanks.txt")
taxonomy <- read_qza("taxonomy.qza")$data %>% parse_taxonomy()
taxonomy$Kingdom <- gsub("d__", "",taxonomy$Kingdom)

metadata <- data.frame(metadata[,-1], row.names=metadata[,1])

taxasums<-summarize_taxa(SVs, taxonomy)$Species

# now I'm going to create betadiversity plots and look at the what environmental
# variables may be influencing beta diversity
mtaxasums <- t(as.matrix(taxasums))
set.seed(123)
nmds <- metaMDS(mtaxasums, distance = "bray")
## Square root transformation
## Wisconsin double standardization
## Run 0 stress 0.1929938 
## Run 1 stress 0.2104758 
## Run 2 stress 0.2122065 
## Run 3 stress 0.2017629 
## Run 4 stress 0.203068 
## Run 5 stress 0.2175578 
## Run 6 stress 0.207389 
## Run 7 stress 0.398229 
## Run 8 stress 0.1902876 
## ... New best solution
## ... Procrustes: rmse 0.07528947  max resid 0.3149976 
## Run 9 stress 0.1891515 
## ... New best solution
## ... Procrustes: rmse 0.03557354  max resid 0.1759746 
## Run 10 stress 0.1891534 
## ... Procrustes: rmse 0.003089313  max resid 0.01691136 
## Run 11 stress 0.1952534 
## Run 12 stress 0.1907372 
## Run 13 stress 0.1952534 
## Run 14 stress 0.1952537 
## Run 15 stress 0.2033445 
## Run 16 stress 0.1888986 
## ... New best solution
## ... Procrustes: rmse 0.01907387  max resid 0.08086317 
## Run 17 stress 0.1909827 
## Run 18 stress 0.1918261 
## Run 19 stress 0.1937858 
## Run 20 stress 0.2114092 
## *** Best solution was not repeated -- monoMDS stopping criteria:
##     20: stress ratio > sratmax
en <- envfit(nmds, metadata, permutations = 999, na.rm = TRUE)

# base r plots are a mess so going to use ggplot2
#extract NMDS scores (x and y coordinates) for sites from newer versions of vegan package
data.scores <- as.data.frame(scores(nmds)$sites)
data.score <- merge(data.scores, metadata, by=0) # adding metadata
#extracting vectors
en_coord_cont = as.data.frame(scores(en, "vectors")) * ordiArrowMul(en)
en_coord_cat = as.data.frame(scores(en, "factors")) * ordiArrowMul(en)

Mantel Tests

# To work out which environmental vectors to include I'm going to do some mantel tests
# first create the distance matrices then do the mantel test
dist.abundance <- vegdist(t(taxasums), method = "bray")

dist.time <- dist(metadata$Time.Since.Start, method = "euclidean")
abund_time = mantel(dist.abundance, dist.time, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_time
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.time, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.5896 
##       Significance: 1e-04 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0600 0.0809 0.1001 0.1265 
## Permutation: free
## Number of permutations: 9999
dist.day.night <- dist(metadata$nDay_Night, method = "euclidean")
abund_day.night = mantel(dist.abundance, dist.day.night, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_day.night
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.day.night, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: -0.06712 
##       Significance: 0.9412 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0605 0.0815 0.1013 0.1267 
## Permutation: free
## Number of permutations: 9999
dist.moon <- dist(metadata$nMoon_Moonset, method = "euclidean")
abund_moon = mantel(dist.abundance, dist.moon, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_moon
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.moon, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: -0.02434 
##       Significance: 0.6761 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0622 0.0850 0.1042 0.1305 
## Permutation: free
## Number of permutations: 9999
dist.tide <- dist(metadata$UKHO_Tide_Height, method = "euclidean")
abund_tide = mantel(dist.abundance, dist.tide, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_tide
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.tide, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.01872 
##       Significance: 0.2957 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0553 0.0755 0.0948 0.1187 
## Permutation: free
## Number of permutations: 9999
dist.caltemp <- dist(metadata$Calista_Temp_oC, method = "euclidean")
abund_caltemp = mantel(dist.abundance, dist.caltemp, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_caltemp
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.caltemp, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.3094 
##       Significance: 1e-04 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0783 0.1022 0.1257 0.1508 
## Permutation: free
## Number of permutations: 9999
dist.airtemp <- dist(metadata$SotonMet_Air_Temp_oC, method = "euclidean")
abund_airtemp = mantel(dist.abundance, dist.airtemp, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_airtemp
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.airtemp, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.02901 
##       Significance: 0.2779 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0694 0.0921 0.1152 0.1457 
## Permutation: free
## Number of permutations: 9999
dist.wind <- dist(metadata$SotonMet_Wind_Speed_knots, method = "euclidean")
abund_wind = mantel(dist.abundance, dist.wind, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_wind
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.wind, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.06955 
##       Significance: 0.1595 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0943 0.1225 0.1476 0.1779 
## Permutation: free
## Number of permutations: 9999
dist.winddirection <- dist(metadata$SotonMet_Wind_Direction_degrees, method = "euclidean")
abund_winddirection = mantel(dist.abundance, dist.winddirection, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_winddirection
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.winddirection, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.4445 
##       Significance: 1e-04 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0944 0.1240 0.1488 0.1825 
## Permutation: free
## Number of permutations: 9999
dist.gust <- dist(metadata$SotonMet_Max_Gust_knots, method = "euclidean")
abund_gust = mantel(dist.abundance, dist.gust, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_gust
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.gust, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.05695 
##       Significance: 0.1869 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0884 0.1174 0.1442 0.1723 
## Permutation: free
## Number of permutations: 9999
dist.pressure <- dist(metadata$SotonMet_Barometric_Pressure_millibars, method = "euclidean")
abund_pressure = mantel(dist.abundance, dist.pressure, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_pressure
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.pressure, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.2172 
##       Significance: 0.0049 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0983 0.1303 0.1606 0.1946 
## Permutation: free
## Number of permutations: 9999
dist.UoStemp <- dist(metadata$UoS_Temp_oC, method = "euclidean")
abund_UoStemp = mantel(dist.abundance, dist.UoStemp, method = "spearman", permutations = 9999, na.rm = TRUE)
abund_UoStemp
## 
## Mantel statistic based on Spearman's rank correlation rho 
## 
## Call:
## mantel(xdis = dist.abundance, ydis = dist.UoStemp, method = "spearman",      permutations = 9999, na.rm = TRUE) 
## 
## Mantel statistic r: 0.3463 
##       Significance: 1e-04 
## 
## Upper quantiles of permutations (null model):
##    90%    95%  97.5%    99% 
## 0.0833 0.1088 0.1306 0.1561 
## Permutation: free
## Number of permutations: 9999

nMDS

#removing non significant variables and time since start
drops <- c("Conc.NOC",
           "Conc.ESS",
           "Conc.PCR2",
           "Conc.diff",
           "nDay_Night",
           "nMoon_Moonset",
           "SampleNo",
           "UKHO_Tide_Height",
           "Calista_Temp_oC",
           "SotonMet_Air_Temp_oC",
           "SotonMet_Wind_Speed_knots",
           "SotonMet_Max_Gust_knots",
           "SotonMet_Wind_Direction_degrees",
           "Time.Since.Start")
en_coord_cont <- en_coord_cont[!(row.names(en_coord_cont) %in% drops),]
en_label=c("P[b]","Temp")

# now to plot and fit the environmental data
nMDS16S = ggplot(data = data.score, aes(x = NMDS1, y = NMDS2)) + 
  geom_point(data = data.score, aes(colour = `Time.Since.Start`/60), size = 3, alpha = 0.8) + 
  scale_fill_viridis(limits = c(0, 94))  +
  scale_colour_viridis(limits = c(0, 94), 
                        breaks = c(0, 24, 48, 72),
                        labels = c(0, 24, 48, 72)) +
  geom_segment(aes(x = 0, y = 0, xend = NMDS1, yend = NMDS2), 
               data = en_coord_cont, linewidth =1, alpha = 0.5, colour = "grey30") +
  geom_text(data = en_coord_cont, aes(x = NMDS1, y = NMDS2), colour = "grey30", 
            fontface = "bold", label = en_label , parse = TRUE, position = position_nudge(x = -0.05)) + 
  theme(axis.title = element_text(size = 10, face = "bold", colour = "grey30"), 
        panel.background = element_blank(), panel.border = element_rect(fill = NA, colour = "grey30"), 
        axis.ticks = element_blank(), axis.text = element_blank(), legend.key = element_blank(), 
        legend.title = element_text(size = 10, face = "bold", colour = "grey30"), 
        legend.text = element_text(size = 9, colour = "grey30")) +
  theme_q2r() +
  labs(title="Prokaryotes", colour = "Hours Since Start")
nMDS16S

## Alpha Diversity Time Series

shannon<-read_qza("/Users/RMS1U18/Library/CloudStorage/OneDrive-UniversityofSouthampton/PhD/Soton/Bioinformatics/16S/core-metrics-results/shannon_vector.qza")
shannon<-shannon$data %>% rownames_to_column("Sample") # this moves the sample names to a new column that matches the metadata and allows them to be merged

observedASVs<-read_qza("/Users/RMS1U18/Library/CloudStorage/OneDrive-UniversityofSouthampton/PhD/Soton/Bioinformatics/16S/core-metrics-results/observed_features_vector.qza")
observedASVs<-observedASVs$data %>% rownames_to_column("Sample")

evenness<-read_qza("/Users/RMS1U18/Library/CloudStorage/OneDrive-UniversityofSouthampton/PhD/Soton/Bioinformatics/16S/core-metrics-results/evenness_vector.qza")
evenness<-evenness$data %>% rownames_to_column("Sample")

faithPD<-read_qza("/Users/RMS1U18/Library/CloudStorage/OneDrive-UniversityofSouthampton/PhD/Soton/Bioinformatics/16S/core-metrics-results/faith_pd_vector.qza")
faithPD<-faithPD$data %>% rownames_to_column("Sample")

#assign alpha diversity metrics to the metadata
metadata<-
  metadata %>% 
  left_join(shannon) %>%
  left_join(observedASVs) %>%
  left_join(evenness) %>%
  left_join(faithPD)
## Joining with `by = join_by(Sample)`
## Joining with `by = join_by(Sample)`
## Joining with `by = join_by(Sample)`
## Joining with `by = join_by(Sample)`
metadata$`Start.Date.Time`<-dmy_hm(metadata$`Start.Date.Time`) # getting the dates into ISO format

ShannonPlot <- ggplot(metadata, aes(x = `Start.Date.Time`, y=`shannon_entropy`)) +
  geom_line() +
  geom_point() +
  ylim(2,7.5) +
  theme_q2r() +
  xlab("Time") +
  ylab("Diversity")
ShannonPlot 

EvennessPlot <- ggplot(metadata, aes(x = `Start.Date.Time`, y=`pielou_evenness`)) +
  geom_line() +
  geom_point() +
  ylim(0.25,0.85) +
  theme_q2r() +
  xlab("Time") +
  ylab("Evenness")
EvennessPlot 

ObservedPlot <- ggplot(metadata, aes(x = `Start.Date.Time`, y=`observed_features`)) +
  geom_line() +
  geom_point() +
  theme_q2r() +
  xlab("Time") +
  ylab("ASVs") +
  labs(title = "Prokaryotes")
ObservedPlot

FaithPlot <- ggplot(metadata, aes(x = `Start.Date.Time`, y=`faith_pd`)) +
  geom_line() +
  geom_point() +
  theme_q2r() +
  xlab("Time") +
  ylab("Faith PD")
FaithPlot 

Alpha16STime<- ggarrange(ObservedPlot, ShannonPlot, EvennessPlot, FaithPlot,
          ncol = 1,
          align = "hv")

Stacked Taxa Time Series

Species

# now I'm going to create a summary of the top 10 taxa and the remainder of taxa
setwd("/Users/RMS1U18/Library/CloudStorage/OneDrive-UniversityofSouthampton/PhD/Soton/Bioinformatics/16S/")
metadata <- read_q2metadata("metadata-no-blanks.txt")
metadata <- data.frame(metadata[,-1], row.names=metadata[,1])
taxasums <- cbind(taxasums, total = rowSums(taxasums)) #add total col for each species
taxasums <- taxasums[order(-taxasums$total),]
write.csv(taxasums, "SpeciesList.csv")
# seperate the top 10 taxa
top10species <- taxasums %>%
  arrange(desc(total)) %>%
  slice(1:10)
#remove the top 10 species 
remainderSpecies <- taxasums[-(1:10),]
remainderSpecies <- rbind(remainderSpecies, Remainder = colSums(remainderSpecies)) #add total row for each sample
# rejoin the tables
taxasums10 <- rbind(top10species, remainderSpecies)
# take the top 11 which will now include the remainder total
taxasums10 <- taxasums10 %>%
  arrange(desc(total)) %>%
  slice(1:11)
#remove the total column
taxasums10$total <- NULL
#rownames to first column
taxasums10<- tibble::rownames_to_column(taxasums10)
taxasums10<- data.frame(t(taxasums10)) #transpose
names(taxasums10)<- taxasums10[1,] # make the taxa the rownames
taxasums10<- taxasums10[-1,] # remove the first duplicate first row

# add Date_Time to taxasums10
t.metadata <- metadata[,c("Sample","Start.Date.Time")]
t.metadata$Start.Date.Time <- as.character(t.metadata$Start.Date.Time) #convert posixct back to character

#t.metadata <- row_to_names(t.metadata, row_number = 1, remove_row = TRUE) #make the sample names col names
taxasums101 <- cbind(t.metadata, taxasums10) #combine taxa and times
taxasums101$Start.Date.Time<- dmy_hm(taxasums101$Start.Date.Time)
write.csv(taxasums101, "species10.csv")

taxasums101$Sample<- NULL # remove sample column

species<- taxasums101 %>%
  pivot_longer(!Start.Date.Time,, names_to = "Taxa", values_to = "Sequence count")


species$`Sequence count`<- as.numeric(species$`Sequence count`)

#plot stack
species16Sstack <- species %>%
  ggplot(aes(fill= Taxa,
             x=Start.Date.Time,
             y=`Sequence count`,
             colour = Taxa)) +
  geom_area(position="stack", 
            stat="identity") +
  geom_line(position = 'stack', linewidth = 1, colour = "black" ) +
  scale_colour_viridis(discrete = TRUE, direction = -1) +
  scale_fill_viridis(discrete = TRUE, direction = -1) +
  theme_q2r() +
  theme(legend.position = "bottom", legend.box = "vertical", legend.justification = "left", legend.box.just = "left", legend.direction = "vertical") +
  labs(title = "Prokaryotes", x="Time")
  

species16Sstack

### Genus

Environmental Time Series

# set directory
setwd("~/OneDrive - University of Southampton/PhD/Soton/Soton_metadata/")
#import UoS temp data
UoStemp2019<-read.csv("UoSbuoy2019Tempdata_Cut.csv",
           stringsAsFactors = FALSE)

#import calista temp data
caltemp<-read.csv("calista_temp_dock_only.csv",
           stringsAsFactors = FALSE)

#import UKHO time an moon data
UKHO<-read.csv("UKHO_tidesunmoondata.csv",
                  stringsAsFactors = FALSE)

#import sotonmet dat inc. tide height (use stringAsFactors= FALSE to stop the data_time being read as a factor)
sotmet11<-read.csv("Sotonmet11Nov2019.csv",
                   stringsAsFactors = FALSE)
sotmet12<-read.csv("Sotonmet12Nov2019.csv",
                   stringsAsFactors = FALSE)
sotmet13<-read.csv("Sotonmet13Nov2019.csv",
                   stringsAsFactors = FALSE)
sotmet14<-read.csv("Sot14Nov2019.csv",
                   stringsAsFactors = FALSE)
sotmet15<-read.csv("Sot15Nov2019.csv",
                   stringsAsFactors = FALSE)

#merge both sotonmet files
sotmet<- rbind(sotmet11,sotmet12,sotmet13,sotmet14,sotmet15)

# convert Date_Time into the same format for each table
UoStemp2019$Date.Time <- dmy_hm(UoStemp2019$Date.Time)
#I'm now replacing the Date.Time column head in UoStemp2019 to Date_Time
names(UoStemp2019)[2]  <- "Date_Time"
UoStemp2019 = subset(UoStemp2019, select = -c(Date,Time) )#xts doesn't work when there are Data Time columns that aren't in the right class

caltemp$Date_Time <- mdy_hms(caltemp$Date_Time)

#in this line I am also pasting the Date and time columns into a single Date_Time coloumn
sotmet$Date_Time <- dmy_hm(paste(sotmet$Date, sotmet$Time))
sotmet = subset(sotmet, select = -c(Date,Time) )

UKHO$Date_Time <- dmy_hm(UKHO$Date_Time)

#Import the data into a time series object

caltemp.ts = as.xts(caltemp, frequency = 60*24) # measurements were taken every minute and I'm setting my season as 1 day
caltemp.ts.interp = na.approx(object = caltemp.ts,
                              xout = seq(min(caltemp$Date_Time), max(caltemp$Date_Time), "min")) # using linear interpolation to fill in the missing values
plot(caltemp.ts.interp$Temperature..Deg.C., col = 2)

lines(caltemp.ts$Temperature..Deg.C., col = 1) #red section shows were the miss data has been filling in

sotmet.ts = as.xts(sotmet, frequency = 12*24) # these observations were every 5 minutes
sotmet.ts.interp = na.approx(object = sotmet.ts,
                              xout = seq(min(sotmet$Date_Time), max(sotmet$Date_Time), "5 min")) #adds approximation to missing values


UKHO.ts = as.xts(UKHO$Tide_Height.m., order.by = UKHO$Date_Time, frequency = 60*24) # these observations were every minute
plot(sotmet.ts.interp$DEPTH, col = 2)

lines(UKHO.ts, col = 1) #comparing tide hieghts data

UoStemp2019.ts = as.xts(UoStemp2019, frequency = 4*24)# these observations were every 15 minutes
plot(UoStemp2019.ts)

lines(caltemp.ts.interp$Temperature..Deg.C., col = 2) #comparing temperature data

# Converting all to 15 min intervals
timeLim<- as.POSIXct(strptime(c("2019-11-11 12:00:00", "2019-11-15 09:00:00"), format = "%Y-%m-%d %H:%M:%S"))
UKHO.15.ts <- rollmean(UKHO.ts, 15, align = "right")
plot(UKHO.ts, col=1)

lines(UKHO.15.ts, col=2) 

UKHO.15.ts<-fortify(UKHO.15.ts)
tidePlot <- ggplot(UKHO.15.ts, aes(x=Index, y=UKHO.15.ts)) +
  geom_line() + 
  scale_x_datetime(limits = timeLim, date_breaks = "1 day", date_labels = "%e %b") +
  xlab("Time") +
  ylab("Tide (m)") +
  theme_q2r()
tidePlot
## Warning: Removed 1605 rows containing missing values (`geom_line()`).

caltemp.i15.ts <- rollmean(caltemp.ts.interp, 15, align = "right")
plot(caltemp.i15.ts$Temperature..Deg.C., col=1)

lines(caltemp.ts.interp$Temperature..Deg.C., col=2)

caltemp.i15.ts <- fortify(caltemp.i15.ts)
tempPlot <- ggplot(caltemp.i15.ts, aes(x=Index, y=Temperature..Deg.C.)) +
  geom_line() + 
  scale_x_datetime(limits = timeLim, date_breaks = "1 day", date_labels = "%e %b") +
  xlab("Time") +
  ylab("Water (°C)") +
  theme_q2r()
tempPlot
## Warning: Removed 1605 rows containing missing values (`geom_line()`).

plot(UoStemp2019.ts, col=2)
lines(caltemp.i15.ts$Temperature..Deg.C., col=1)

sotATMP.15.ts <- rollmean(sotmet.ts.interp$ATMP, 3, align = "right")
plot(sotmet.ts.interp$ATMP, col=1)

lines(sotATMP.15.ts, col=2)

sotATMP.15.ts <- fortify(sotATMP.15.ts)
airTempPlot <- ggplot(sotATMP.15.ts, aes(x=Index, y=ATMP)) +
  geom_line() + 
  scale_x_datetime(limits = timeLim, date_breaks = "1 day", date_labels = "%e %b") +
  xlab("Time") +
  ylab("Air (°C)") +
  theme_q2r()
airTempPlot
## Warning: Removed 321 rows containing missing values (`geom_line()`).

sotWSPD.15.ts <- rollmean(sotmet.ts.interp$WSPD, 3, align = "right")
plot(sotmet.ts.interp$WSPD, col=1)

lines(sotWSPD.15.ts, col=2)

sotWSPD.15.ts <- fortify(sotWSPD.15.ts)
windPlot <- ggplot(sotWSPD.15.ts, aes(x=Index, y=WSPD)) +
  geom_line() + 
  scale_x_datetime(limits = timeLim, date_breaks = "1 day", date_labels = "%e %b") +
  xlab("Time") +
  ylab("Wind (knots)") +
  theme_q2r()
windPlot
## Warning: Removed 321 rows containing missing values (`geom_line()`).

sotWD.15.ts <- rollmean(sotmet.ts.interp$WD, 3, align = "right")
plot(sotmet.ts.interp$WD, col=1)

lines(sotWD.15.ts, col=2)

sotGST.15.ts <- rollmean(sotmet.ts.interp$GST, 3, align = "right")
plot(sotmet.ts.interp$GST, col=1)

lines(sotGST.15.ts, col=2)

gustPlot <- ggplot(sotGST.15.ts, aes(x=Index, y=GST)) +
  geom_line() + 
  scale_x_datetime(limits = timeLim, date_breaks = "1 day", date_labels = "%e %b") +
  xlab("Time") +
  ylab("Gust (knots)") +
  theme_q2r()
gustPlot
## Warning: Removed 321 rows containing missing values (`geom_line()`).

sotBARO.15.ts <- rollmean(sotmet.ts.interp$BARO, 3, align = "right")
plot(sotmet.ts.interp$BARO, col=1)

lines(sotBARO.15.ts, col=2)

sotBARO.15.ts <- fortify(sotBARO.15.ts)
presurePlot <- ggplot(sotBARO.15.ts, aes(x=Index, y=BARO)) +
  geom_line() + 
  scale_x_datetime(limits = timeLim, date_breaks = "1 day", date_labels = "%e %b") +
  xlab("Time") +
  ylab(expression(P[b]*" (mbars)")) +
  theme_q2r()
presurePlot
## Warning: Removed 321 rows containing missing values (`geom_line()`).

# Plot together

envPlot <- ggarrange(tidePlot, tempPlot, airTempPlot, presurePlot, windPlot,
                     ncol = 1,
                     align = "hv")
## Warning: Removed 1605 rows containing missing values (`geom_line()`).
## Removed 1605 rows containing missing values (`geom_line()`).
## Warning: Removed 321 rows containing missing values (`geom_line()`).
## Removed 321 rows containing missing values (`geom_line()`).
## Removed 321 rows containing missing values (`geom_line()`).
envPlot

ggsave("EnvTimeSeries.pdf", height = 5, width = 3.5, device = "pdf")

Combined Plots

Beta Diversity

nMDSPlot<- ggarrange(nMDS16S, nMDS18S,
                      nrow = 1,
                      align = "hv",
                      common.legend = TRUE,
                      legend = "bottom")
nMDSPlot

ggsave("nMDS.pdf", height=4, width=7, device="pdf")

Alpha Diversity Time Series

AlphaTime<- ggarrange(Alpha16STime, Alpha18STime,
                      ncol = 2,
                      align = "hv")
AlphaTime

ggsave("AlphaDiversity.pdf", height=4.7, width=7, device="pdf")

Stacked Taxa Time Series

speciesStacks<- ggarrange(species16Sstack, species18SStack, 
                          ncol = 1,
                          align = "hv", 
                          legend = "bottom")
speciesStacks

ggsave("Specie10TimeSeries.pdf", height=10, width=6, device="pdf")

Tables

library("sjPlot")
## Install package "strengejacke" from GitHub (`devtools::install_github("strengejacke/strengejacke")`) to load all sj-packages at once!
setwd("~/OneDrive - University of Southampton/PhD/Soton/Bioinformatics/16S")
alphaMeta<- read_tsv("alpha16S-metadata.tsv")
## Rows: 39 Columns: 24
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: "\t"
## dbl (24): Sample, Time-Since-Start, UKHO_Tide_Height, Calista_Temp_oC, Soton...
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
alphaEnvCorr<- tab_corr(alphaMeta, corr.method = "spearman", show.p = TRUE, p.numeric = TRUE)
alphaEnvCorr
  Sample Time-Since-Start UKHO_Tide_Height Calista_Temp_oC SotonMet_Air_Temp_oC SotonMet_Wind_Speed_knots SotonMet_Wind_Direction_degrees SotonMet_Max_Gust_knots SotonMet_Barometric_Pressure_millibars UoS_Temp_oC faith_pd pielou_evenness observed_features shannon_entropy Bacteria; Proteobacteria;
Alphaproteobacteria; SAR11_clade;
Clade_I; Clade_Ia; NA
Bacteria; Proteobacteria;
Alphaproteobacteria; Rhodobacterales;
Rhodobacteraceae; Planktomarina; NA
Bacteria; Bacteroidota; Bacteroidia;
Flavobacteriales; Cryomorphaceae;
uncultured; NA
Bacteria; Proteobacteria;
Gammaproteobacteria; Pseudomonadales;
Thioglobaceae; SUP05_cluster; NA
Bacteria; Bacteroidota; Bacteroidia;
Flavobacteriales; Flavobacteriaceae;
NS5_marine_group; NA
Bacteria; Proteobacteria;
Gammaproteobacteria; Pseudomonadales;
Saccharospirillaceae; Oleispira; NA
Bacteria; Proteobacteria;
Gammaproteobacteria; Burkholderiales;
Methylophilaceae; OM43_clade; NA
Bacteria; Bacteroidota; Bacteroidia;
Flavobacteriales; Flavobacteriaceae;
NS3a_marine_group; uncultured_bacterium
Bacteria; Proteobacteria;
Gammaproteobacteria; Enterobacterales;
Alteromonadaceae; Glaciecola; NA
Bacteria; Proteobacteria;
Alphaproteobacteria; Puniceispirillales;
SAR116_clade;
Candidatus_Puniceispirillum;
uncultured_bacterium
Sample   -1.000
(<.001)
0.000
(1.00)
0.805
(<.001)
0.257
(.115)
0.198
(.227)
0.436
(.006)
0.207
(.207)
0.076
(.645)
0.883
(<.001)
0.025
(.879)
0.274
(.092)
0.075
(.649)
0.233
(.152)
0.512
(.001)
0.004
(.979)
0.657
(<.001)
0.411
(.009)
0.306
(.058)
-0.695
(<.001)
0.362
(.024)
0.655
(<.001)
0.496
(.001)
0.532
(<.001)
Time-Since-Start -1.000
(<.001)
  0.000
(1.00)
-0.805
(<.001)
-0.257
(.115)
-0.198
(.227)
-0.436
(.006)
-0.207
(.207)
-0.076
(.645)
-0.883
(<.001)
-0.025
(.879)
-0.274
(.092)
-0.075
(.649)
-0.233
(.152)
-0.512
(.001)
-0.004
(.979)
-0.657
(<.001)
-0.411
(.009)
-0.306
(.058)
0.695
(<.001)
-0.362
(.024)
-0.655
(<.001)
-0.496
(.001)
-0.532
(<.001)
UKHO_Tide_Height 0.000
(1.00)
0.000
(1.00)
  0.204
(.212)
0.024
(.887)
0.106
(.520)
0.228
(.162)
0.052
(.754)
0.021
(.901)
-0.198
(.227)
-0.150
(.361)
-0.221
(.177)
-0.072
(.662)
-0.219
(.180)
-0.018
(.911)
-0.109
(.509)
-0.290
(.073)
0.226
(.167)
0.062
(.709)
0.031
(.853)
0.184
(.262)
-0.300
(.064)
-0.161
(.326)
0.112
(.496)
Calista_Temp_oC 0.805
(<.001)
-0.805
(<.001)
0.204
(.212)
  0.426
(.007)
0.116
(.482)
0.487
(.002)
0.118
(.473)
0.057
(.730)
0.727
(<.001)
-0.155
(.343)
0.280
(.085)
-0.132
(.423)
0.008
(.963)
0.371
(.021)
-0.013
(.939)
0.456
(.004)
0.244
(.135)
0.187
(.254)
-0.486
(.002)
0.177
(.281)
0.472
(.003)
0.272
(.094)
0.338
(.035)
SotonMet_Air_Temp_oC 0.257
(.115)
-0.257
(.115)
0.024
(.887)
0.426
(.007)
  0.368
(.021)
0.227
(.165)
0.386
(.015)
-0.216
(.187)
0.442
(.005)
0.015
(.926)
0.262
(.108)
-0.172
(.295)
0.023
(.888)
0.223
(.172)
-0.041
(.805)
0.091
(.582)
0.059
(.720)
-0.129
(.435)
-0.125
(.447)
-0.027
(.872)
0.010
(.953)
-0.081
(.624)
-0.146
(.376)
SotonMet_Wind_Speed_knots 0.198
(.227)
-0.198
(.227)
0.106
(.520)
0.116
(.482)
0.368
(.021)
  -0.022
(.895)
0.961
(<.001)
0.217
(.185)
0.260
(.110)
0.094
(.570)
-0.352
(.028)
0.017
(.917)
-0.132
(.423)
-0.069
(.674)
0.367
(.021)
0.253
(.120)
0.386
(.015)
-0.130
(.429)
-0.348
(.030)
-0.054
(.743)
0.120
(.468)
0.279
(.085)
0.131
(.428)
SotonMet_Wind_Direction_degrees 0.436
(.006)
-0.436
(.006)
0.228
(.162)
0.487
(.002)
0.227
(.165)
-0.022
(.895)
  0.030
(.857)
0.060
(.717)
0.296
(.067)
0.080
(.627)
-0.002
(.989)
-0.025
(.881)
0.019
(.911)
0.432
(.006)
-0.049
(.769)
0.217
(.186)
0.260
(.110)
0.305
(.059)
-0.471
(.002)
0.292
(.071)
0.157
(.341)
0.033
(.841)
0.347
(.031)
SotonMet_Max_Gust_knots 0.207
(.207)
-0.207
(.207)
0.052
(.754)
0.118
(.473)
0.386
(.015)
0.961
(<.001)
0.030
(.857)
  0.218
(.182)
0.292
(.071)
0.100
(.546)
-0.318
(.049)
-0.045
(.785)
-0.152
(.355)
-0.036
(.829)
0.389
(.014)
0.294
(.069)
0.372
(.020)
-0.064
(.701)
-0.367
(.022)
-0.062
(.708)
0.191
(.244)
0.359
(.025)
0.130
(.429)
SotonMet_Barometric_Pressure_millibars 0.076
(.645)
-0.076
(.645)
0.021
(.901)
0.057
(.730)
-0.216
(.187)
0.217
(.185)
0.060
(.717)
0.218
(.182)
  -0.060
(.715)
-0.459
(.003)
-0.182
(.267)
-0.318
(.048)
-0.357
(.026)
0.054
(.742)
0.456
(.003)
0.282
(.082)
0.437
(.005)
0.201
(.220)
-0.495
(.001)
0.038
(.818)
0.150
(.363)
0.333
(.038)
0.274
(.091)
UoS_Temp_oC 0.883
(<.001)
-0.883
(<.001)
-0.198
(.227)
0.727
(<.001)
0.442
(.005)
0.260
(.110)
0.296
(.067)
0.292
(.071)
-0.060
(.715)
  0.089
(.589)
0.301
(.063)
0.018
(.915)
0.195
(.235)
0.389
(.014)
0.064
(.697)
0.598
(<.001)
0.231
(.157)
0.124
(.450)
-0.582
(<.001)
0.148
(.367)
0.586
(<.001)
0.441
(.005)
0.323
(.045)
faith_pd 0.025
(.879)
-0.025
(.879)
-0.150
(.361)
-0.155
(.343)
0.015
(.926)
0.094
(.570)
0.080
(.627)
0.100
(.546)
-0.459
(.003)
0.089
(.589)
  -0.125
(.447)
0.572
(<.001)
0.474
(.003)
-0.032
(.845)
-0.121
(.461)
-0.101
(.541)
-0.235
(.150)
-0.072
(.662)
0.217
(.185)
-0.040
(.806)
0.017
(.916)
-0.125
(.448)
0.043
(.795)
pielou_evenness 0.274
(.092)
-0.274
(.092)
-0.221
(.177)
0.280
(.085)
0.262
(.108)
-0.352
(.028)
-0.002
(.989)
-0.318
(.049)
-0.182
(.267)
0.301
(.063)
-0.125
(.447)
  -0.166
(.312)
0.377
(.019)
0.375
(.019)
-0.229
(.161)
0.178
(.277)
-0.040
(.811)
0.174
(.288)
-0.023
(.890)
0.238
(.144)
0.277
(.088)
0.051
(.756)
0.132
(.424)
observed_features 0.075
(.649)
-0.075
(.649)
-0.072
(.662)
-0.132
(.423)
-0.172
(.295)
0.017
(.917)
-0.025
(.881)
-0.045
(.785)
-0.318
(.048)
0.018
(.915)
0.572
(<.001)
-0.166
(.312)
  0.804
(<.001)
-0.019
(.910)
-0.126
(.444)
-0.023
(.890)
0.082
(.618)
-0.082
(.621)
0.059
(.723)
-0.045
(.787)
-0.022
(.895)
-0.018
(.913)
0.196
(.232)
shannon_entropy 0.233
(.152)
-0.233
(.152)
-0.219
(.180)
0.008
(.963)
0.023
(.888)
-0.132
(.423)
0.019
(.911)
-0.152
(.355)
-0.357
(.026)
0.195
(.235)
0.474
(.003)
0.377
(.019)
0.804
(<.001)
  0.198
(.226)
-0.186
(.256)
0.135
(.410)
0.132
(.423)
0.021
(.901)
-0.022
(.893)
0.119
(.471)
0.152
(.355)
0.059
(.723)
0.273
(.093)
Bacteria; Proteobacteria;
Alphaproteobacteria; SAR11_clade;
Clade_I; Clade_Ia; NA
0.512
(.001)
-0.512
(.001)
-0.018
(.911)
0.371
(.021)
0.223
(.172)
-0.069
(.674)
0.432
(.006)
-0.036
(.829)
0.054
(.742)
0.389
(.014)
-0.032
(.845)
0.375
(.019)
-0.019
(.910)
0.198
(.226)
  -0.493
(.001)
0.095
(.564)
0.426
(.007)
0.693
(<.001)
-0.390
(.014)
0.820
(<.001)
0.146
(.374)
0.017
(.919)
0.481
(.002)
Bacteria; Proteobacteria;
Alphaproteobacteria; Rhodobacterales;
Rhodobacteraceae; Planktomarina; NA
0.004
(.979)
-0.004
(.979)
-0.109
(.509)
-0.013
(.939)
-0.041
(.805)
0.367
(.021)
-0.049
(.769)
0.389
(.014)
0.456
(.003)
0.064
(.697)
-0.121
(.461)
-0.229
(.161)
-0.126
(.444)
-0.186
(.256)
-0.493
(.001)
  0.351
(.028)
0.144
(.381)
-0.290
(.073)
-0.317
(.049)
-0.517
(.001)
0.214
(.191)
0.332
(.039)
0.044
(.792)
Bacteria; Bacteroidota; Bacteroidia;
Flavobacteriales; Cryomorphaceae;
uncultured; NA
0.657
(<.001)
-0.657
(<.001)
-0.290
(.073)
0.456
(.004)
0.091
(.582)
0.253
(.120)
0.217
(.186)
0.294
(.069)
0.282
(.082)
0.598
(<.001)
-0.101
(.541)
0.178
(.277)
-0.023
(.890)
0.135
(.410)
0.095
(.564)
0.351
(.028)
  0.396
(.013)
0.056
(.735)
-0.755
(<.001)
-0.007
(.967)
0.903
(<.001)
0.817
(<.001)
0.392
(.014)
Bacteria; Proteobacteria;
Gammaproteobacteria; Pseudomonadales;
Thioglobaceae; SUP05_cluster; NA
0.411
(.009)
-0.411
(.009)
0.226
(.167)
0.244
(.135)
0.059
(.720)
0.386
(.015)
0.260
(.110)
0.372
(.020)
0.437
(.005)
0.231
(.157)
-0.235
(.150)
-0.040
(.811)
0.082
(.618)
0.132
(.423)
0.426
(.007)
0.144
(.381)
0.396
(.013)
  0.284
(.080)
-0.689
(<.001)
0.497
(.001)
0.228
(.162)
0.468
(.003)
0.685
(<.001)
Bacteria; Bacteroidota; Bacteroidia;
Flavobacteriales; Flavobacteriaceae;
NS5_marine_group; NA
0.306
(.058)
-0.306
(.058)
0.062
(.709)
0.187
(.254)
-0.129
(.435)
-0.130
(.429)
0.305
(.059)
-0.064
(.701)
0.201
(.220)
0.124
(.450)
-0.072
(.662)
0.174
(.288)
-0.082
(.621)
0.021
(.901)
0.693
(<.001)
-0.290
(.073)
0.056
(.735)
0.284
(.080)
  -0.348
(.030)
0.740
(<.001)
0.243
(.137)
0.170
(.301)
0.565
(<.001)
Bacteria; Proteobacteria;
Gammaproteobacteria; Pseudomonadales;
Saccharospirillaceae; Oleispira; NA
-0.695
(<.001)
0.695
(<.001)
0.031
(.853)
-0.486
(.002)
-0.125
(.447)
-0.348
(.030)
-0.471
(.002)
-0.367
(.022)
-0.495
(.001)
-0.582
(<.001)
0.217
(.185)
-0.023
(.890)
0.059
(.723)
-0.022
(.893)
-0.390
(.014)
-0.317
(.049)
-0.755
(<.001)
-0.689
(<.001)
-0.348
(.030)
  -0.318
(.049)
-0.604
(<.001)
-0.667
(<.001)
-0.536
(<.001)
Bacteria; Proteobacteria;
Gammaproteobacteria; Burkholderiales;
Methylophilaceae; OM43_clade; NA
0.362
(.024)
-0.362
(.024)
0.184
(.262)
0.177
(.281)
-0.027
(.872)
-0.054
(.743)
0.292
(.071)
-0.062
(.708)
0.038
(.818)
0.148
(.367)
-0.040
(.806)
0.238
(.144)
-0.045
(.787)
0.119
(.471)
0.820
(<.001)
-0.517
(.001)
-0.007
(.967)
0.497
(.001)
0.740
(<.001)
-0.318
(.049)
  0.045
(.783)
-0.001
(.993)
0.494
(.001)
Bacteria; Bacteroidota; Bacteroidia;
Flavobacteriales; Flavobacteriaceae;
NS3a_marine_group; uncultured_bacterium
0.655
(<.001)
-0.655
(<.001)
-0.300
(.064)
0.472
(.003)
0.010
(.953)
0.120
(.468)
0.157
(.341)
0.191
(.244)
0.150
(.363)
0.586
(<.001)
0.017
(.916)
0.277
(.088)
-0.022
(.895)
0.152
(.355)
0.146
(.374)
0.214
(.191)
0.903
(<.001)
0.228
(.162)
0.243
(.137)
-0.604
(<.001)
0.045
(.783)
  0.793
(<.001)
0.457
(.003)
Bacteria; Proteobacteria;
Gammaproteobacteria; Enterobacterales;
Alteromonadaceae; Glaciecola; NA
0.496
(.001)
-0.496
(.001)
-0.161
(.326)
0.272
(.094)
-0.081
(.624)
0.279
(.085)
0.033
(.841)
0.359
(.025)
0.333
(.038)
0.441
(.005)
-0.125
(.448)
0.051
(.756)
-0.018
(.913)
0.059
(.723)
0.017
(.919)
0.332
(.039)
0.817
(<.001)
0.468
(.003)
0.170
(.301)
-0.667
(<.001)
-0.001
(.993)
0.793
(<.001)
  0.438
(.005)
Bacteria; Proteobacteria;
Alphaproteobacteria; Puniceispirillales;
SAR116_clade;
Candidatus_Puniceispirillum;
uncultured_bacterium
0.532
(<.001)
-0.532
(<.001)
0.112
(.496)
0.338
(.035)
-0.146
(.376)
0.131
(.428)
0.347
(.031)
0.130
(.429)
0.274
(.091)
0.323
(.045)
0.043
(.795)
0.132
(.424)
0.196
(.232)
0.273
(.093)
0.481
(.002)
0.044
(.792)
0.392
(.014)
0.685
(<.001)
0.565
(<.001)
-0.536
(<.001)
0.494
(.001)
0.457
(.003)
0.438
(.005)
 
Computed correlation used spearman-method with listwise-deletion.
setwd("~/OneDrive - University of Southampton/PhD/Soton/Bioinformatics/18S_dada2redo")
alphaMeta<- read_tsv("alpha18S-metadata.tsv")
## Rows: 39 Columns: 24
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: "\t"
## dbl (24): id, Time-Since-Start, UKHO_Tide_Height, Calista_Temp_oC, SotonMet_...
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
alphaEnvCorr<- tab_corr(alphaMeta, corr.method = "spearman", show.p = TRUE, p.numeric = TRUE)
alphaEnvCorr
  id Time-Since-Start UKHO_Tide_Height Calista_Temp_oC SotonMet_Air_Temp_oC SotonMet_Wind_Speed_knots SotonMet_Wind_Direction_degrees SotonMet_Max_Gust_knots SotonMet_Barometric_Pressure_millibars UoS_Temp_oC faith_pd pielou_evenness observed_features shannon_entropy Eukaryota; NA; NA; NA; NA; NA; NA Eukaryota; Arthropoda; Maxillopoda; NA;
NA; NA; NA
Eukaryota; Annelida; Polychaeta;
Capitellida; Capitellida; Capitellida;
NA
Eukaryota; Dinoflagellata; Dinophyceae;
NA; NA; NA; NA
Eukaryota; Cryptophyceae; Cryptophyceae;
Cryptomonadales; Cryptomonadales;
uncultured; Cryptophyta_sp.
Eukaryota; Picozoa; Picomonadida;
Picomonadida; Picomonadida;
Picomonadida; uncultured_phototrophic
Eukaryota; Diatomea; Mediophyceae;
Mediophyceae; Mediophyceae;
Thalassiosira; NA
Eukaryota; Protalveolata; Syndiniales;
Syndiniales; Syndiniales_Group_I;
Syndiniales_Group_I;
Karlodinium_veneficum
Eukaryota; Arthropoda; Maxillopoda;
Calanoida; Calanoida; Calanoida;
uncultured_eukaryote
Eukaryota; Incertae_Sedis;
Incertae_Sedis; Incertae_Sedis;
Incertae_Sedis; Telonema;
Telonema_antarcticum
id   -1.000
(<.001)
0.000
(1.00)
0.805
(<.001)
0.257
(.115)
0.198
(.227)
0.436
(.006)
0.207
(.207)
0.076
(.645)
0.883
(<.001)
-0.048
(.772)
-0.305
(.060)
-0.349
(.030)
-0.278
(.087)
0.050
(.764)
-0.168
(.305)
0.314
(.052)
0.195
(.233)
-0.309
(.056)
0.028
(.864)
-0.209
(.201)
-0.225
(.168)
0.339
(.035)
-0.031
(.852)
Time-Since-Start -1.000
(<.001)
  0.000
(1.00)
-0.805
(<.001)
-0.257
(.115)
-0.198
(.227)
-0.436
(.006)
-0.207
(.207)
-0.076
(.645)
-0.883
(<.001)
0.048
(.772)
0.305
(.060)
0.349
(.030)
0.278
(.087)
-0.050
(.764)
0.168
(.305)
-0.314
(.052)
-0.195
(.233)
0.309
(.056)
-0.028
(.864)
0.209
(.201)
0.225
(.168)
-0.339
(.035)
0.031
(.852)
UKHO_Tide_Height 0.000
(1.00)
0.000
(1.00)
  0.204
(.212)
0.024
(.887)
0.106
(.520)
0.228
(.162)
0.052
(.754)
0.021
(.901)
-0.198
(.227)
-0.001
(.995)
0.046
(.783)
-0.079
(.635)
-0.020
(.903)
0.140
(.394)
0.185
(.259)
-0.078
(.638)
-0.088
(.593)
0.025
(.880)
0.069
(.677)
0.215
(.190)
0.113
(.491)
0.372
(.020)
0.238
(.144)
Calista_Temp_oC 0.805
(<.001)
-0.805
(<.001)
0.204
(.212)
  0.426
(.007)
0.116
(.482)
0.487
(.002)
0.118
(.473)
0.057
(.730)
0.727
(<.001)
-0.325
(.044)
-0.246
(.131)
-0.544
(<.001)
-0.327
(.043)
0.159
(.331)
-0.155
(.343)
0.353
(.028)
0.115
(.483)
-0.455
(.004)
0.111
(.500)
-0.121
(.462)
-0.008
(.962)
0.520
(.001)
0.053
(.749)
SotonMet_Air_Temp_oC 0.257
(.115)
-0.257
(.115)
0.024
(.887)
0.426
(.007)
  0.368
(.021)
0.227
(.165)
0.386
(.015)
-0.216
(.187)
0.442
(.005)
-0.151
(.360)
-0.151
(.357)
-0.361
(.024)
-0.249
(.127)
0.118
(.474)
0.046
(.782)
0.133
(.420)
0.080
(.628)
-0.288
(.075)
0.030
(.857)
-0.093
(.572)
0.170
(.301)
0.262
(.107)
-0.039
(.814)
SotonMet_Wind_Speed_knots 0.198
(.227)
-0.198
(.227)
0.106
(.520)
0.116
(.482)
0.368
(.021)
  -0.022
(.895)
0.961
(<.001)
0.217
(.185)
0.260
(.110)
0.097
(.557)
-0.202
(.218)
-0.166
(.312)
-0.220
(.179)
0.021
(.897)
-0.001
(.997)
0.028
(.866)
0.078
(.638)
-0.043
(.796)
0.036
(.828)
-0.107
(.518)
-0.220
(.178)
-0.082
(.622)
-0.144
(.383)
SotonMet_Wind_Direction_degrees 0.436
(.006)
-0.436
(.006)
0.228
(.162)
0.487
(.002)
0.227
(.165)
-0.022
(.895)
  0.030
(.857)
0.060
(.717)
0.296
(.067)
-0.122
(.459)
-0.012
(.941)
-0.298
(.065)
-0.113
(.492)
0.041
(.804)
0.042
(.799)
0.027
(.869)
0.403
(.011)
-0.014
(.934)
-0.031
(.850)
-0.190
(.247)
-0.079
(.632)
0.475
(.002)
-0.039
(.815)
SotonMet_Max_Gust_knots 0.207
(.207)
-0.207
(.207)
0.052
(.754)
0.118
(.473)
0.386
(.015)
0.961
(<.001)
0.030
(.857)
  0.218
(.182)
0.292
(.071)
0.177
(.280)
-0.224
(.171)
-0.198
(.227)
-0.237
(.146)
0.042
(.800)
-0.014
(.930)
0.079
(.632)
0.106
(.520)
-0.054
(.746)
0.045
(.784)
-0.110
(.504)
-0.187
(.255)
-0.045
(.787)
-0.120
(.467)
SotonMet_Barometric_Pressure_millibars 0.076
(.645)
-0.076
(.645)
0.021
(.901)
0.057
(.730)
-0.216
(.187)
0.217
(.185)
0.060
(.717)
0.218
(.182)
  -0.060
(.715)
-0.232
(.155)
0.157
(.341)
-0.220
(.178)
0.108
(.513)
-0.247
(.130)
0.033
(.842)
-0.121
(.464)
0.294
(.069)
0.041
(.803)
-0.095
(.564)
-0.132
(.424)
-0.313
(.053)
0.079
(.632)
-0.105
(.525)
UoS_Temp_oC 0.883
(<.001)
-0.883
(<.001)
-0.198
(.227)
0.727
(<.001)
0.442
(.005)
0.260
(.110)
0.296
(.067)
0.292
(.071)
-0.060
(.715)
  -0.078
(.638)
-0.319
(.047)
-0.402
(.011)
-0.312
(.053)
-0.013
(.937)
-0.094
(.571)
0.310
(.055)
0.049
(.767)
-0.361
(.024)
-0.093
(.575)
-0.209
(.202)
-0.256
(.116)
0.283
(.080)
-0.160
(.330)
faith_pd -0.048
(.772)
0.048
(.772)
-0.001
(.995)
-0.325
(.044)
-0.151
(.360)
0.097
(.557)
-0.122
(.459)
0.177
(.280)
-0.232
(.155)
-0.078
(.638)
  0.202
(.216)
0.678
(<.001)
0.386
(.016)
-0.082
(.618)
-0.086
(.600)
-0.110
(.502)
0.081
(.624)
0.385
(.016)
-0.050
(.761)
0.027
(.871)
-0.054
(.744)
-0.288
(.076)
0.101
(.541)
pielou_evenness -0.305
(.060)
0.305
(.060)
0.046
(.783)
-0.246
(.131)
-0.151
(.357)
-0.202
(.218)
-0.012
(.941)
-0.224
(.171)
0.157
(.341)
-0.319
(.047)
0.202
(.216)
  0.381
(.017)
0.934
(<.001)
-0.073
(.659)
-0.011
(.949)
-0.198
(.227)
0.171
(.296)
0.391
(.014)
-0.090
(.585)
0.050
(.764)
0.038
(.816)
-0.080
(.629)
0.100
(.544)
observed_features -0.349
(.030)
0.349
(.030)
-0.079
(.635)
-0.544
(<.001)
-0.361
(.024)
-0.166
(.312)
-0.298
(.065)
-0.198
(.227)
-0.220
(.178)
-0.402
(.011)
0.678
(<.001)
0.381
(.017)
  0.634
(<.001)
-0.087
(.598)
-0.091
(.583)
-0.341
(.034)
-0.066
(.688)
0.452
(.004)
-0.118
(.476)
0.000
(1.00)
-0.043
(.797)
-0.473
(.002)
-0.009
(.957)
shannon_entropy -0.278
(.087)
0.278
(.087)
-0.020
(.903)
-0.327
(.043)
-0.249
(.127)
-0.220
(.179)
-0.113
(.492)
-0.237
(.146)
0.108
(.513)
-0.312
(.053)
0.386
(.016)
0.934
(<.001)
0.634
(<.001)
  -0.092
(.577)
-0.067
(.684)
-0.209
(.201)
0.125
(.448)
0.438
(.005)
-0.105
(.525)
0.038
(.819)
0.002
(.990)
-0.221
(.176)
0.085
(.609)
Eukaryota; NA; NA; NA; NA; NA; NA 0.050
(.764)
-0.050
(.764)
0.140
(.394)
0.159
(.331)
0.118
(.474)
0.021
(.897)
0.041
(.804)
0.042
(.800)
-0.247
(.130)
-0.013
(.937)
-0.082
(.618)
-0.073
(.659)
-0.087
(.598)
-0.092
(.577)
  -0.456
(.004)
0.016
(.925)
0.056
(.733)
-0.427
(.007)
0.846
(<.001)
0.455
(.004)
0.789
(<.001)
0.233
(.153)
0.668
(<.001)
Eukaryota; Arthropoda; Maxillopoda; NA;
NA; NA; NA
-0.168
(.305)
0.168
(.305)
0.185
(.259)
-0.155
(.343)
0.046
(.782)
-0.001
(.997)
0.042
(.799)
-0.014
(.930)
0.033
(.842)
-0.094
(.571)
-0.086
(.600)
-0.011
(.949)
-0.091
(.583)
-0.067
(.684)
-0.456
(.004)
  -0.178
(.277)
-0.095
(.563)
0.011
(.945)
-0.522
(.001)
-0.199
(.224)
-0.323
(.045)
-0.150
(.361)
-0.226
(.167)
Eukaryota; Annelida; Polychaeta;
Capitellida; Capitellida; Capitellida;
NA
0.314
(.052)
-0.314
(.052)
-0.078
(.638)
0.353
(.028)
0.133
(.420)
0.028
(.866)
0.027
(.869)
0.079
(.632)
-0.121
(.464)
0.310
(.055)
-0.110
(.502)
-0.198
(.227)
-0.341
(.034)
-0.209
(.201)
0.016
(.925)
-0.178
(.277)
  0.016
(.925)
-0.040
(.808)
0.123
(.456)
0.012
(.940)
0.011
(.947)
0.149
(.366)
-0.088
(.594)
Eukaryota; Dinoflagellata; Dinophyceae;
NA; NA; NA; NA
0.195
(.233)
-0.195
(.233)
-0.088
(.593)
0.115
(.483)
0.080
(.628)
0.078
(.638)
0.403
(.011)
0.106
(.520)
0.294
(.069)
0.049
(.767)
0.081
(.624)
0.171
(.296)
-0.066
(.688)
0.125
(.448)
0.056
(.733)
-0.095
(.563)
0.016
(.925)
  0.340
(.034)
0.069
(.676)
-0.229
(.162)
0.017
(.918)
0.003
(.986)
0.140
(.394)
Eukaryota; Cryptophyceae; Cryptophyceae;
Cryptomonadales; Cryptomonadales;
uncultured; Cryptophyta_sp.
-0.309
(.056)
0.309
(.056)
0.025
(.880)
-0.455
(.004)
-0.288
(.075)
-0.043
(.796)
-0.014
(.934)
-0.054
(.746)
0.041
(.803)
-0.361
(.024)
0.385
(.016)
0.391
(.014)
0.452
(.004)
0.438
(.005)
-0.427
(.007)
0.011
(.945)
-0.040
(.808)
0.340
(.034)
  -0.398
(.012)
-0.049
(.768)
-0.353
(.028)
-0.386
(.015)
-0.217
(.185)
Eukaryota; Picozoa; Picomonadida;
Picomonadida; Picomonadida;
Picomonadida; uncultured_phototrophic
0.028
(.864)
-0.028
(.864)
0.069
(.677)
0.111
(.500)
0.030
(.857)
0.036
(.828)
-0.031
(.850)
0.045
(.784)
-0.095
(.564)
-0.093
(.575)
-0.050
(.761)
-0.090
(.585)
-0.118
(.476)
-0.105
(.525)
0.846
(<.001)
-0.522
(.001)
0.123
(.456)
0.069
(.676)
-0.398
(.012)
  0.488
(.002)
0.810
(<.001)
0.110
(.505)
0.778
(<.001)
Eukaryota; Diatomea; Mediophyceae;
Mediophyceae; Mediophyceae;
Thalassiosira; NA
-0.209
(.201)
0.209
(.201)
0.215
(.190)
-0.121
(.462)
-0.093
(.572)
-0.107
(.518)
-0.190
(.247)
-0.110
(.504)
-0.132
(.424)
-0.209
(.202)
0.027
(.871)
0.050
(.764)
0.000
(1.00)
0.038
(.819)
0.455
(.004)
-0.199
(.224)
0.012
(.940)
-0.229
(.162)
-0.049
(.768)
0.488
(.002)
  0.539
(<.001)
-0.118
(.476)
0.438
(.005)
Eukaryota; Protalveolata; Syndiniales;
Syndiniales; Syndiniales_Group_I;
Syndiniales_Group_I;
Karlodinium_veneficum
-0.225
(.168)
0.225
(.168)
0.113
(.491)
-0.008
(.962)
0.170
(.301)
-0.220
(.178)
-0.079
(.632)
-0.187
(.255)
-0.313
(.053)
-0.256
(.116)
-0.054
(.744)
0.038
(.816)
-0.043
(.797)
0.002
(.990)
0.789
(<.001)
-0.323
(.045)
0.011
(.947)
0.017
(.918)
-0.353
(.028)
0.810
(<.001)
0.539
(<.001)
  0.117
(.477)
0.756
(<.001)
Eukaryota; Arthropoda; Maxillopoda;
Calanoida; Calanoida; Calanoida;
uncultured_eukaryote
0.339
(.035)
-0.339
(.035)
0.372
(.020)
0.520
(.001)
0.262
(.107)
-0.082
(.622)
0.475
(.002)
-0.045
(.787)
0.079
(.632)
0.283
(.080)
-0.288
(.076)
-0.080
(.629)
-0.473
(.002)
-0.221
(.176)
0.233
(.153)
-0.150
(.361)
0.149
(.366)
0.003
(.986)
-0.386
(.015)
0.110
(.505)
-0.118
(.476)
0.117
(.477)
  0.009
(.958)
Eukaryota; Incertae_Sedis;
Incertae_Sedis; Incertae_Sedis;
Incertae_Sedis; Telonema;
Telonema_antarcticum
-0.031
(.852)
0.031
(.852)
0.238
(.144)
0.053
(.749)
-0.039
(.814)
-0.144
(.383)
-0.039
(.815)
-0.120
(.467)
-0.105
(.525)
-0.160
(.330)
0.101
(.541)
0.100
(.544)
-0.009
(.957)
0.085
(.609)
0.668
(<.001)
-0.226
(.167)
-0.088
(.594)
0.140
(.394)
-0.217
(.185)
0.778
(<.001)
0.438
(.005)
0.756
(<.001)
0.009
(.958)
 
Computed correlation used spearman-method with listwise-deletion.